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Abstract 

In recent years extensive studies on the Earth’s climate system have been 
carried out by means of advanced complex network statistics. The great 
majority of these studies, however, have been focusing on investigating cor¬ 
relation structures within single climatic fields directly on or parallel to the 
Earth’s surface. Here, we develop a novel approach of node weighted coupled 
network measures to study correlations between ocean and atmosphere in the 
Northern Hemisphere extratropics and construct 18 coupled climate networks, 
each consisting of two subnetworks. In all cases, one subnetwork represents 
monthly sea-surface temperature (SST) anomalies, while the other is based on 
the monthly geopotential height (HGT) of isobaric surfaces at different pres¬ 
sure levels covering the troposphere as well as the lower stratosphere. The 
weighted cross-degree density proves to be consistent with the leading coupled 
pattern obtained from maximum covariance analysis. Network measures of 
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higher order allow for a further analysis of the correlation structure between 
the two fields and consistently indicate that in the Northern Hemisphere ex¬ 
tratropics the ocean is correlated with the atmosphere in a hierarchical fashion 
such that large areas of the ocean surface correlate with multiple statistically 
dissimilar regions in the atmosphere. Ultimately we show that, this observed 
hierarchy is linked to large-scale atmospheric variability patterns, such as the 
Pacific North American pattern, forcing the ocean on monthly time scales. 

Keywords: coupled climate networks, extratropical ocean-atmosphere inter¬ 
action, node-weighted network measures, hierarchical networks 


1 Introduction 


In the last years, complex network analysis has been established as a powerful tool 


to study statistical interdependencies in the climate system (Donges et al. 2009a 


Donges et al. 2015aJ Tsonis and Roebber 2004; Tsonis et al. 2006; Tsonis et al. 
2008). Links in the so-called climate networks represent functional interdependencies 


indicated by significant correlations (Donges et ah 2009a; Donges et ah 2009b; Palus 


et ah 2011; Radebach et ah 2013) or the synchronous occurrence of extreme events 
(Boers et ah 2013; Boers et ah 2014b; Malik et ah 2010; Malik et ah 2011; Stolbova 


et ah 2014) in climatic time series taken at different grid points or measurement 


sites on or parallel to the Earth’s surface. 

In addition to studies on observational data of climate dynamics, climate net¬ 
works have also been applied successfully to hindcast extreme events, such as ex¬ 
treme precipitation in South America (Boers et ah 2014a), or to predict the oc¬ 
currence of El Nino episodes (Ludescher et ah 2013; Ludescher et ah 2014) and 


discriminate between different event types (Radebach et ah 2013; Tsonis and Swan¬ 


son 


2008; Wiedermann et ah 2016; Yamasaki et ah 2008). 


So far, most studies conducted within the framework of climate networks focused 
solely on the dynamics within a single climatic field or layer. Besides atmospheric 
characteristics like surface air temperature, sea level pressure, or precipitation, re¬ 
cent studies have also addressed ocean dynamics represented by ocean temperature 


variability at the surface (Feng and Dijkstra 2014; Tantet and Dijkstra 2014) or 
different depths (Mheen et ah 2013) as well as the spatio-temporal variability in the 
strength of the Atlantic meridional overturning circulation (Feng et ah 2014). 

It is well known, however, that the dynamics within the two major subcompo¬ 
nents of the Earth’s climate system, ocean and atmosphere, are closely entangled 


(Frankignoul et ah 2001; Trenberth and Hurrell 1994). Examples for these inter¬ 
relationships include the North Atlantic eddy-driven jet stream (Woollings et ah 


2010) or the Pacific ocean forcing to the atmosphere which is closely related to the 


dynamics of the El Nino Southern Oscillation (Wyrtki 1975). Further, it has been 
shown that on time scales of up to one month the ocean is forced by atmospheric 
circulation, prominently manifested in terms of long-term variability patterns like 


the Pacific North American pattern (e.g. Frankignoul and Sennchael 2007) and the 


North Atlantic Oscillation (Czaja and Frankignoul 1999; Gastineau and Frankignoul 


2015). 
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Inspired by approaches to investigate the interaction structure between different 


mutually coupled subsystems such as infrastructure networks (Boccaletti et ah 2014 


Buldyrev et al. 2010; Vespignani 2010) a novel set of coupled network measures has 
been proposed by Donges et al. (2011) which provides a general tool to quantify 
interdependencies between subcomponents in complex coupled climate networks. 
The latter framework has been successfully applied to investigate interactions be¬ 
tween different layers of geopotential height fields, where each isobaric surface forms 
a subcomponent of a larger climate network. Similarly, coupled climate networks 
have been constructed to study ocean-atmosphere interactions in the tropical Pa¬ 
cific (Feng et al. 2012) or over the South Atlantic Convergence Zone (Tirabassi et al. 


2015). 


Following upon these previous studies, in this work we extend the approach 
by Donges et al. (2011) and present an exploratory study to understand and quantify 
ocean-atmosphere interactions in the Northern Hemisphere mid-to-high latitudes 
during boreal winter at monthly scales. This temporal restriction is chosen, since 
previous studies by means of lagged maximum covariance analysis (MCA) have 
already revealed that the statistical interrelationship between atmosphere and ocean 
is strongest and most significant at lags of zero or one month during late fall and 
winter (e.g. Czaja and Frankignoul 1999; Frankignoul and Sennchael 2007; Gastineau 


and Frankignoul 2015J; Liu et al. 2006; Wen et al. 2005). 

To investigate further the spatial structure of these complex interaction patterns, 
we construct here in total 18 coupled climate networks consisting of two layers 
each, one layer representing sea surface temperature (SST) anomalies and the other 
geopotential height fields (HGT) at different pressure levels from 1000 to 10 mbar 
covering the entire troposphere as well as the lower stratosphere. 

Our area of study covers the whole Northern Hemisphere north of 30°N so that 
the density of grid points in the considered climate data sets increases rapidly to¬ 
wards the poles and induces some bias in the unweighted network measures (Rade- 
bach et al. 2013; Tsonis et al. 2006). Therefore, the standard coupled network ap¬ 
proach by Donges et al. (2011) is not sufficient in the present case. To overcome the 
problem associated with the heterogeneous spatial density of grid points interpreted 
as nodes of the climate network, Heitzig et al. (2012) introduced a novel set of net¬ 
work measures that takes into account the different sizes or weights of nodes in the 
network. By following an axiomatic approach, each standard (unweighted) network 
measure can be transformed into its weighted counterpart, the so-called node split¬ 
ting invariant (n.s.i.) network measure. Corresponding n.s.i. measures have also 
been derived by Zemp et al. (2014) for edge-weighted and directed networks. 

To quantify the topology of coupled climate networks, we rely in this work on the 
previously defined versions of local (i.e. node-wise) n.s.i. coupled network measures 
(Feng et al. 2012; Wiedermann et al. 2013) and additionally derive further weighted 
global network measures following the approach introduced by Heitzig et al. (2012). 
This allows us to assess and compare the macroscopic correlation structure in each 
of the 18 coupled climate networks. 

We compare the results of MCA (e.g. Storch and Zwiers[2001), a well-established 
standard tool from statistical climatology, with the cross-degree density of nodes in 
the different subnetworks and confirm expected similarities between the two mea- 
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sures (Donges et al. 2015a). By utilizing network measures of higher order such as 
the n.s.i. local cross-clustering coefficient, we find that the statistical interrelation 
between ocean and atmosphere exhibits a hierarchical structure, in which individual 
parts or areas of the ocean surface correlate strongly with multiple statistically dis¬ 
similar parts of the atmosphere. Building upon previous studies by, e.g., Czaja and 
Frankignoul 1999; Frankignoul and Sennchael 2007, and Gastineau and Frankignoul 
2015 we relate the observed hierarchy to dominant atmospheric patterns forcing the 
ocean on the time scales investigated in this study. 

The remainder of this paper is organized as follows. Section [2] introduces the 
data sets and all methods, i.e., maximum covariance analysis and coupled climate 
network analysis, that are applied in this study. Section [3] presents all results of the 
analysis followed by conclusions and an outlook discussing future research tasks in 
Section 01 


2 Data & Methods 


2.1 Data description 

We construct coupled climate networks from two different climatic observables in 
order to investigate their interaction structure. One subnetwork is based on monthly 
anomalies of geopotential height (HGT) fields obtained from the ERA40 reanalysis 
project of the European Centre for Medium-Range Weather Forecast (Uppala et 


al. 2005). The data is given on a regular latitude/longitude grid with a spatial 
resolution of AA = A <fi = 2.5°. In total, we investigate 18 layers of HGT fields. 
The corresponding pressure values at each isobaric surface as well as the average 
geopotential height are given in Tab. [l] The second subnetwork is constructed from 
the monthly averaged SST field (HadlSSTl) provided by the Met Office Hadley 
Centre (Rayner et al. 2003) with a resolution of AA = A0 = 1°. All grid points 
with corresponding time series containing missing values are removed from the data 
set as they represent areas which have been at least temporarily covered by sea-ice. 

For our analysis we investigate all grid points north of A = 30° N excluding the 
North Pole itself. Both data sets are cropped in their temporal extent to cover 
the same time span from January 1958 to December 2001 and, hence, each time 
series consists of T = 528 temporal sampling points. We obtain a total number 
of N s = 6201 grid points for the SST data and W = 3456 grid points for each 
isobaric surface i of HGT. For both data sets, we remove the annual cycle by 
subtracting the climatic mean for each month from each time series. Since we 
focus on the spatial structure of strong statistical interrelationships between ocean 
and atmosphere during boreal winter months (DJF), we use only the corresponding 
values which yields a length of each time series of r = 132 data points. 


2.2 Maximum covariance analysis (MCA) 

Consider two sets of time series {X Sn (t)}^ : 1 1 and {Aj m (f)}^ =1 representing two dif¬ 
ferent climatic fields, which in the scope of our application are the SST field (in what 
follows indicated by the index s ) and one layer i of HGT (see also Tab. [I]). Further 
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assume each individual time series in both fields to be normalized to zero mean and 
unit variance. The linear lag-zero cross-correlation matrix C r( with entries C Snhn is 
then defined as 


= ( 1 ) 

T ti 

where n = 1,... ,N S , m = 1,... ,Ni and r denotes the total number of temporal 
sampling points in the two time series. Due to the heterogeneous spatial distribution 
of grid points in the present data sets all matrix entries C Sn i m are additionally 
multiplied by the square roots of the cosine of latitudinal positions A. to ensure 
equal weighting. This then yields the weighted cross-correlation matrix C"- with 
entries 


c Zim = \/cos A Sn cos X^C Snim . 


( 2 ) 


Analogously to empirical orthogonal function (EOF) analysis (e.g. Ghil et al. 2002; 


Hannachi et al. 


(p. 


(m) 

si 


2007), MCA identifies orthonormal pairs of coupled patterns /y m 1 = 

m) 




■ Ps$) and = {Pif 1 ■ ■ .p y C.) f° r m — 1 ,... ,R (with R being the rank 
of C sl ) which explain as much as possible of the covariance between pairs of time 
series taken from the two different climatic fields (e.g. Bretherton et al. 1992; Storch 
and Zwiers 2001). The coupled patterns are obtained by solving the singular value 
problem of the weighted cross-covariance matrix, 


c "pf” 1 = (3) 

(C» T p < /” ) = aJP- (4) 

They are ordered according to their respective singular values cq with cq > cr 2 > 
... > < 7 f { . Hence, cq denotes the largest among the R singular values that can 
be found to solve the above equations. Therefore, p: i l> and f>s ) are referred to as 
the leading coupled patterns representing the largest fraction of squared covariance 
between the two climatic fields given by 

2.3 Coupled climate network construction 

In climate networks, each node represents a climatic time series and links indicate 
strong correlations between two series. Hence, the N x N (N = N s + Nj) correlation 
matrix contains the pairwise linear statistical relationships between all time series 
considered for the network construction. Here, we independently construct coupled 
climate networks for all combinations of the SST field and each of the 18 isobaric sur¬ 
faces of HGT, which shall be investigated separately and rely on the linear Pearson 
correlation coefficient as an appropriate measure of statistical association. Hence, 
the correlation matrix has the form 
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c = 




(5) 


The two block matrices C ss (N s x N s ) and C vl (N t x N t ) represent the (internal) 
correlation matrices of the SST and HGT fields, respectively, which consist of ele¬ 
ments 


l T 

Cs n Sm = ~ ^ ] X Sn (t) X Srrl (t) , Tl, 171 = 1 , . . . , N s , 

T t! 

i t 

C inim = n,m = l,...,Ni. 

T ti 


( 6 ) 

(7) 


The elements of C st = Cf s are derived according to Eq. ([Tj) . Note, that (in contrast to 
the computation of the leading coupled patterns) we construct the coupled climate 
networks from the unweighted correlation matrix C, while the correction for the 
heterogeneous spatial distribution of nodes is implemented into the corresponding 
network measures (see Sec. 2.4). 


From the correlation matrix C one generally derives the network’s adjacency ma¬ 
trix A + by setting a fixed threshold T such that only a certain fraction (i.e. the fink 
density p) of strongest correlations is represented by finks in the resulting climate 
network. For obtaining the adjacency matrix A + of coupled climate networks, we 
refine this procedure by fixing a desired link density p s = pi = 0.01 for the structure 
of internal links within the two subnetworks representing SST and HGT fields, re¬ 
spectively. This means that only nodes with a correlation above the empirical 99th 
percentile of correlations between all time series within each field are connected. 
This condition then leads to internal correlation thresholds T s = 0.8101 for the SST 
field and Tj for each isobaric surface of HGT (Tab. [I]). Usually, the dynamics within 
the different climatic fields shows much higher cross-correlations than between them. 
We account for this fact by assuming the fraction of significant interactions between 
the climatic fields to be lower than those within them. Specifically, we request a 
cross-link density of p sl = 0.005 < p s = Pi, which is lower than the internal ones, and 
derive a cross-threshold T s j for each layer of HGT individually (Fig. [I]). All internal 
thresholds T s and Tj are significantly larger than the obtained cross-thresholds T s j. 
Thus, setting a global fink density or threshold would cause no or few cross-links 
to be present between the two fields or respective subnetworks. We further note 
that all links in each of the coupled climate networks represent correlations that 
are significant at least at the 95% confidence level of a standard t-test, where the 
degrees of freedom are determined by the total number of temporal sampling points 
r in each of the time series (with r = 132 we thus obtain 130 degrees of freedom 
when neglecting the presence of serial correlations in the individual time series). 

The different values of T si already give an impression of the strength of correla¬ 
tion between the SST field and the different isobaric layers: low thresholds generally 
indicate weaker correlations while high thresholds imply stronger similarity between 
both fields. Further we note that the resulting cross-thresholds vary smoothly with 
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the choice of cross-link density (Fig. [Tj) and we thus consider the construction mech¬ 
anism to be sufficiently insensitive to the actual choice of cross-link density. 

Using the different thresholds introduced above, we obtain the coupled climate 
network’s adjacency matrix by individually thresholding the absolute correlation 
values between and within both fields as 


f 0(|C SS | — T s ) ©(|C si | — T si )\ 
\Q(\C is \-T si ) ©(|C Ui\-TJJ 


where ©(•) denotes the Heaviside function. Note that in most recent studies on 
climate networks self-loops (resulting in a non-vanishing trace of the adjacency ma¬ 
trix) have been excluded. In this case the adjacency matrix is usually denoted as 
A. Since we aim to apply node splitting invariant network measures (see below) to 
quantify the network’s topology we specifically demand each node to be connected 
with itself. The resulting matrix A + is referred to as the extended adjacency ma¬ 
trix (Heitzig et al. 2012). Further note, that the usage of the term coupled in coupled 
climate networks does not imply the notion of any directionality or causal influence 
between the two fields under study. It is simply meant to indicate the fact that the 
network under study is composed of more than a single climatic field. 


2.4 Coupled network characteristics 

The local (point-wise) and global structure of a climate network can be quantified 


by a variety of network measures (Albert and Barabasi 2002; Donges et al. 2009a 


Newman 2003), which can generally be interpreted as specific operations on the 


adjacency matrix. The climate networks in this study are constructed from climate 
data sets where the density of grid points and, hence, the density of nodes in the 
network, rapidly increases towards the North pole. In order to avoid a bias in the 
evaluation of the climate network’s structure, we account for this effect by relying 
on node-weighted network measures and value nodes with a gradually decreasing 
weight as one moves from the equator towards the pole. To quantify the correlation 
structure between ocean and atmosphere at each node we focus on two previously de¬ 


fined node weighted local network measures, the n.s.i cross-degree (Feng et al. 2012) 
and the n.s.i. local cross-clustering coefficient (Wiedermann et al. 2013). In addi¬ 
tion, we utilize the construction mechanism introduced by Heitzig et al. (2012) to 


convert global interacting network measures (Donges et al. 2011) into their weighted 
counterparts. 


2.4.1 Preliminaries 

Consider a coupled climate network G = (V, E) with a set of nodes V, links E and 
the number of nodes N = \V\. Following the general naming convention in the 
climate network framework we identify every node v £ V with a natural number 
p = 1,..., N, such that p serves as the label of the node as well as an index to cor¬ 
responding network characteristics. The network G is represented by its adjacency 
matrix A with A pq = 1 if (p, q ) £ E, A pq = 0 if (p, q ) ^ E. In this study, each cou¬ 
pled climate network is composed of two subnetworks, G s = (V S ,E SS ) representing 
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the ocean and Gi = (V*, E u ) representing a specific atmospheric layer. The set of 
nodes V divides into subsets V s and V, such that each node belongs to exactly one 
subnetwork (i.e. V = V s U V t and Vgfl V; = 0). Likewise, the set of links E splits into 
internal link sets E ss and E n (connecting nodes within a subnetwork) and cross-link 
sets E si connecting nodes v £ V s with nodes q E Vi in the subnetworks G s and G x , 
respectively (Donges et al. 2011). 

In the present case (as for all regularly gridded climate data sets) the share on 
the entire area of the surface that is represented by each node is governed by its 
latitudinal position A„ on the grid. Following Tsonis et al. (2006), we therefore 
assign to each node v in the climate network a weight 


w v = cos A„. 


( 8 ) 


Note that, by following this convention the climate networks’ node weights w v exhibit 
the same dimension as the weights of the cross-covariance matrix in Eq. Q. 

Heitzig et al. (2012) introduced novel set of node splitting invariant (n.s.i.) 
network measures to quantify the topology of a climate network with such a het¬ 
erogeneous spatial node density for the case of a single-layer network and, hence, 
only one climate variable under study. In fact, the n.s.i. network measures are not 
restricted to climate networks but can be utilized to study any type of single-layer 
complex network where nodes represent entities of different weights. Heitzig et al. 
(2012) further showed that each complex network measure can be transformed into 
its weighted counterpart by using a four-step construction mechanism: 


(a) Sum up weights w v whenever the unweighted measure counts nodes. 

(b) Treat every node v £ V as connected with itself. 

(c) Allow equality in summations over indices v and q wherever the original mea¬ 
sure involves a sum over distinct nodes v and q. 

(d) “Plug in” n.s.i. versions of measures wherever they are used in the definition 
of other measures. 


From the definition of the adjacency matrix A + in Eq. (|8j) we note that step (b) of the 
above scheme is in our case already fulfilled. Wiedermann et al. (2013) and Zemp et 
al. (2014) recently utilized the proposed scheme to convert local interacting network 
measures as well as measures for directed networks into their weighted counterparts. 
Here, we additionally derive n.s.i. versions of some global cross-network measures 
that were introduced by Donges et al. (2011). 


2.4.2 Local measures 

For quantifying local cross-network interactions in coupled climate networks we rely 
on two measures, n.s.i. cross-degree k J v * and n.s.i. local cross-clustering coefficient 
Cl *, that were introduced by Wiedermann et al. (2013) and (for the case of the n.s.i. 
cross-degree) by Feng et al. (2012). These two measures are defined as 











(9) 


H* = ^2 W <l A vqi 
qev j 

Cl * = -h E K A m A > q ™ P e [o, 1]. (10) 

\kv ) p,qeVj 

In contrast to the unweighted cross-degree 


K = E A 

q&Vj 


( 11 ) 


which simply counts nodes q E Vj that are connected with v E V, k 3 * is proportional 
to the share on the considered overall ice-free ocean or isobaric surface area, respec¬ 
tively, that is connected with nodes v E Vj in the other subnetwork. It therefore 
gives a notion of how similar the dynamics at a node v E Vi is to that of the other 
climate variable observed at all available grid points. 

Similar to k 3 ,*, C 3 V * no longer relies on the counting of distinct fully connected 
node triples in the network (as for the classical local clustering coefficient (Newman 


the two subnetworks. It gives the probability that an area represented by a node 
v E Vi is connected with two mutually connected and, hence, dynamically similar, 
areas in the opposite subnetwork. In this spirit, C 3 * estimates how likely areas in 
two different climatic fields or subsystems form clusters of statistical equivalence 
between them. A local accumulation of such connected triples represents clusters of 
closely connected nodes and, in the spirit of climate networks, strongly correlated 
regions. 

In order to make the n.s.i. cross-degree k 3 * comparable between the two subnet¬ 
works, we normalize it by the maximum possible weight that nodes v E V) can be 
connected with, 


2003])) but on the weighted sum of occurrences of triples of connected areas within 


k: 


j* _ 


^qeVj W qAt, 


vq 


W, 


E [0, 1]. 


( 12 ) 


In the spirit of earlier work by Donges et al. (2012) and Donner et al. (2010), 


q ev 3 


w n 


we refer to this quantity as the n.s.i. cross-degree density. Here, Wj = Yl 
denotes the total weight of all nodes q E Vj. For the case of a single-layer network, 
a measure similar to the n.s.i. cross-degree density has been introduced by Tsonis 


et al. (2006) in terms of the area weighted connectivity , which quantities the share 


on the subdomain of interest represented by the entire network G that is connected 
with any nodes v E V. 


Generally, Wiedermann et al. (2013) and Zemp et al. (2014) showed that the 


weighted local cross-network measures improve the representation of a network’s 
topology with inhomogeneous node density within the domain of interest in com¬ 
parison with its unweighted counterparts. 
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Donges et al. (2015a) showed that for the unweighted case cross-degree and 


leading coupled patterns display strong similarity if the first coupled patterns explain 
a high fraction of the system’s covariance. A similar assessment can be made for 
the similarity between the leading coupled patterns obtained from a weighted cross¬ 
covariance matrix and the n.s.i. cross-degree (see supporting information). 


2.4.3 Global measures 

In addition to local (per node) network measures we also aim to characterize the 
macroscopic interaction structure of each pair of coupled climate networks by means 
of global network properties. For coupled climate networks a variety of unweighted 


measures have been proposed by Donges et al. (2011). Here, we utilize the construe 


tion mechanism by Heitzig et al. (2012) to convert two of them into their weighted 
counterparts as well. 


N.s.i. global cross-clustering coefficient. The global cross-clustering coeffi¬ 
cient Cij of a subnetwork Gi gives the probability that for a randomly chosen node 
v G Vi one finds neighbors p,q £ Vj that are mutually linked. It is defined as the 
arithmetic mean of all local cross-clustering coefficients C J V , 


e ‘> = w. E c i- < 13 > 

* v&Vi 

This measure can be converted into its n.s.i. counterpart by calculating the weighted 
mean of all values of C{,*, 


= (Mi 

1 veVi 

Again, analogously to the interpretation of the local n.s.i. measures, C*j no longer 
only measures pure node-wise triangular structures but takes into account the share 
on the Earth’s surface areas involved in the formation of triangular structures. Gen¬ 
erally, large values of C* :I (which are induced by a dominance of connected triples 
between the two subnetworks under consideration) indicate strong transitivity in 
the underlying correlation structure. 

N.s.i. cross-transitivity. The cross-transitivity 7p gives the probability that two 
randomly drawn nodes p,q £ Vj are connected if they have a common neighbor 
v £ Vi. It is given as 


T = 

' K) 


E 


V AAA 

veVi Z^p^qeVj ^vp^pq^qv 


E 


veVi 'Yhp^qE.Vj AypAqv 


(15) 


Like Cij , the cross-transitivity is a measure of organization with respect to the cross¬ 


correlation structure in a coupled network (Donges et al. 2011). However, in contrast 
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to Cij , Tij takes into account the increasing influence of nodes with high cross-degree 
and weighs them more heavily than nodes with low cross-degree. More importantly 
it ignores nodes with no links into the opposite layer, since these nodes display a 
zero cross-degree. The node-weighted variant of 7p can be written as 


_ 

'ij ~ 


E 


£ Pl9 ev$ w v A+ p w p A+ q w q A± 


qv 


E 


veVi E P)9 eVj w vA+ p w P iv q At, 


qv 


EueV; w u(^'v ) Q> 

E^^ff ' 


(16) 


We note that both C*j and similarly evaluate the transitivity of correlations be¬ 
tween the two climatic variables under study and, hence, quantify a similar network 
property. They are derived, however, in a disjoint manner. One measure, C* ;] is 
computed as the weighted average taken over C 3 V *. In contrast, despite suggestions 


by Radebach et al. (2013) to decompose the global transitivity into local contribu¬ 


tions, the n.s.i. cross-transitivity 7~* is defined solely as a global network measure 
with no direct local counterpart. It is important to note that n.s.i. cross-transitivity 
and n.s.i. global cross-clustering coefficient are commonly asymmetric in the sense 
that T*j % Tj* and C* 3 % C*j . 

3 Results 


3.1 Maximum covariance analysis (MCA) 

We start our analysis by computing the leading coupled patterns between the SST 
field and the 18 HGT layers for boreal winter (DJF). Figure [2] displays the results 
for three representative layers of HGT at 50 mbar, 100 mbar and 500 mbar. 

By applying MCA, we detect coherent large-scale patterns of winter SST, which 
co-vary with the winter atmospheric circulation structures instantaneously. The 
leading MCA patterns explain rather large amounts of 42%, 63% and 70% (for the 
500, 100 and 50 mbar pressure level, respectively) of the squared covariance. At 
all levels, the leading MCA mode displays significant SST anomalies over the North 
Pacific with maximum values along the sub-Arctic front near 40°N, and anomalies of 
the opposite sign along the western coast of North America (Fig. [|A,C,E) (An and 
Wang 2005; Frankignoul and Sennchael 2007). Over the North Atlantic, a dipole 


structure is seen between the northern part of the Gulf Stream and the Atlantic 
Ocean south of Greenland including parts of the Davis Strait and the North Atlantic 
current. This pattern resembles the first SST EOF for the Northern Hemisphere 
during boreal winter (not shown). 

This general SST pattern is co-varying with a pressure anomaly pattern showing 
a hemispheric annular-like structure in the upper troposphere and lower stratosphere 
(Fig. [2^3, D). In the mid-troposphere (Fig. HI) , this pattern displays wave-like devia¬ 
tions from the annular structure, which show distinct similarities with the wave-train 
structure of the Pacific North American (PNA) pattern. Therefore, the leading MCA 
mode relates negative SST anomalies along the sub-Arctic front with a positive PNA 
phase. 
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The second MCA mode (not shown, explaining 13%, 17% and 21% of the squared 
covariance fraction for the 500, 100 and 50 mbar pressure level, respectively) dis¬ 
plays the strongest SST anomalies over the North Atlantic. Over that region, the 
SST pattern resembles the northern part of the North Atlantic SST tripole pat¬ 
tern which is related with the North Atlantic Oscillation (NAO) (e.g. Czaja and 


Frankignoul 1999; Czaja and Frankignoul 2002; Gastineau and Frankignoul 2015). 


Accordingly, the co-varying atmospheric pattern in the middle troposphere shows 
the cold ocean/warm land (COWL) pattern (introduced by Wallace et al. (1996)) 
including a NAO-like dipole over the North Atlantic. At higher levels, the co-varying 
atmospheric patterns display a pronounced wave number-2 pattern. 

By applying lagged MCA between SST and mid-tropospheric circulation fields, 
several studies for the North Atlantic and the North Pacific have shown that the 
squared covariance fraction is strongest and most significant at lags of 0 and 1 


month during late fall and winter (e.g. Czaja and Frankignoul 1999; Frankignoul 
and Sennchael 2007; Gastineau and Frankignoul 2015} Liu et al. 2006; Wen et al. 
2005). This points to the forcing of the SST by the dominant atmospheric pattern, 


which is the PNA pattern over the Pacific-North American sector (e.g. Frankignoul 


and Sennchael 2007) and the North Atlantic Oscillation (NAO) over the North 


Atlantic-European region (Czaja and Frankignoul 1999; Gastineau and Frankignoul 


2015). On the other hand, results of lagged MCA analyses with the ocean leading by 


1 to 4 months in Frankignoul and Sennchael (2007) and Gastineau and Frankignoul 


(2015) suggest that the SST anomalies have a substantial influence on the large-scale 


atmospheric circulation at these time-scales. 


3.2 Local coupled network measures 

In order to first demonstrate the general consistency of coupled climate network anal¬ 
ysis in comparison with MCA, we continue by generating coupled climate networks 
between the SST field and the three previously considered layers of geopotential 
height (500 mbar, 100 mbar, 50 mbar). The n.s.i. cross-degree densities k 1 * and 
k** are expected to display similar spatial structures as the corresponding leading 
coupled patterns (Donges et al. 2015a| since the latter explain a high share of the 
cross-covariance between both fields (see supporting information). 

As demonstrated in Fig. [3j the results for «** and k 1 * indeed match well the 
results obtained from the MCA when comparing the locations of maximum values in 
the coupled network’s n.s.i. cross-degree densities to those of maximum or minimum 
values in the leading mode of the MCA. Note, that the n.s.i. cross-degree densities 
and k 1 * take, per definition, only positive values, while coupled patterns display 
both, positive and negative values. Hence, k s * and n s * only reproduce structures 
that coincide with the absolute values of the leading coupled patterns. However, as 
only a certain percentage of squared covariance is explained by the leading coupled 
patterns, we also note differences between the patterns revealed by the two methods. 
In particular, the negative center of action around the North Pole that is detected 
by MCA is only weakly present in the cross-degree density fields k s * for the 50 and 
100 mbar HGT fields (compare Fig. [2)3, D with Fig.[3}3,D). For the ocean, preferably 
the marked structures in the leading coupled patterns in both the Atlantic and the 
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Pacific are well recovered by the cross-degree density n l * while some of the weaker 
structures, e.g. in the Black Sea, are missing. 

Network analysis, however, allows us to undertake a further in-depth analysis of 
the correlation structure between the different layers beyond the information pro¬ 
vided by MCA. The n.s.i. local cross-clustering coefficients C l * and C** (Eq. (®» 
give the probabilities that the dynamics at a grid point in, e.g., the SST field is 
similar to that at two grid points in the HGT field, which behave themselves sta¬ 
tistically similar. Note that in the scope of this work we do not account for any 
possible effects induced by a common external forcing of the fields, which might ar¬ 
tificially induce correlations and, hence, cause the presence of spurious links between 
nodes or triples of nodes. We also do not account for indirect (partial) correlations 
or common driver effects within each of the fields when constructing the coupled 
climate networks. Conditioning out these possible influences by means of informa¬ 
tion theoretic approaches (Runge et al. 2012; Runge et al. 2014) and causal effect 
networks (Kretschmer et al. 2016; Runge et al. 2015) thus remains as a subject of 
future research. 

Figure [I] presents the results for the n.s.i. local cross-clustering coefficients C** 
for nodes in the SST field (Fig. |4]A,C,E) and C s * for nodes in the HGT fields 
(Fig. [4j3,D,F). Most nodes in the SST field tend to display a low n.s.i. local cross¬ 
clustering coefficient C 1 * < 0.2 (Fig# ,C,E) and, thus, preferentially correlate with 
nodes in the HGT fields that are mutually dissimilar and therefore disconnected 
(Fig# In contrast, many nodes in the HGT fields exhibit a comparatively high or 
intermediate n.s.i. local cross-clustering coefficient 0.4 < C®* < 1 (for one of the most 
prominent examples compare nodes located at or above the Pacific in Fig, [ffjB.D) . 
Quantitatively, for the combination of the SST and the 500 mbar HGT field we find 
an n.s.i. global cross-clustering coefficient (Eq. (0) of C* si = 0.16 for SST nodes 
and C* s = 0.28 for 500 mbar HGT nodes. Ignoring those nodes in the averaging that 
display zero n.s.i. cross-degree density we obtain values of C* : = 0.42 and C* s = 0.52 
(note that this definition is different from the one presented in Eq. (14) as we specif¬ 
ically exclude the contribution of nodes with no links to the opposite subnetwork). 
The n.s.i. cross-transitivity (Eq. ([16])) which weighs nodes according to their n.s.i. 
cross-degree density gives values of 7j* = 0.2 and T* s = 0.25 for nodes in ocean 
and atmosphere, respectively. For all three measures, the values computed for the 
atmospheric subnetwork exceed those for the ocean and thus consistently imply that 
nodes in the ocean are less likely to connect with mutually connected nodes in the 
atmosphere than vice versa. 

To further quantify the asymmetries in the correlation structure between ocean 
and atmosphere, we investigate for each node with a given n.s.i. cross-degree den¬ 
sity its corresponding n.s.i. local cross-clustering coefficient in a coupled climate 
network composed of the SST and 500 mbar HGT fields (Fig. [6]). This layer is 
chosen as it provides a good indication of the atmospheric circulation over the area 
of study (Gastineau and Frankignoul 2015; Kushnir et al. 2002). Furthermore, it 


displays among the highest values of T si according to Fig. [TJ which has similarly 
been described as a strong statistical signal by Frankignoul and Sennchael (2007). 

For nodes in the SST field (Fig. [6]A.), we find that C**(k**) tends to follow a 
power-law, C l * ~ (k**) “, which indicates a hierarchical network structure (Ravasz 
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et al. 2002; Ravasz and Barabasi 2003) which, in contrast, is absent for nodes in 


the HGT field (Fig. [6)3). Here, the term hierarchical implies that nodes in the SST 
field strongly correlate with disconnected clusters of statistically similar nodes in the 
HGT field as depicted in Fig. [5| This deduction is further supported by the fact that 
for the HGT field, the distribution of combinations of C°* and k** is more widely 
spread and Cl* generally takes higher values than C l *. This implies that nodes in 
the HGT field show a stronger tendency to correlate with mutually connected nodes 
in the SST field, which can be assumed to display a strong statistical similarity 
among themselves (Molkenthin et al. 2014; Tupikina et al. 2014). To test for the 
robustness of our results we have carried out the same analysis as presented in 
Fig.© for internal link densities of p ss = pa = 0.02 and p ss = p tl = 0.05, and 
corresponding cross-link densities p s i = 0.01 and p m = 0.025 (see Figs. 1 and 2 
in the supporting information). Even though the power-law exponent a slightly 
decreases towards zero with increasing link densities, we find that the qualitative 
findings remain unchanged. We thus consider our analysis to be sufficiently robust 
with respect to the actual choice of link densities. 

As a remark, we note a general tendency of nodes at the boundaries of a cluster 
that links with the opposite field to display comparatively low values of n.s.i. cross¬ 
degree density and increased values of n.s.i. local cross-clustering coefficient (Fig. [3] 
and Fig. [4]) . In contrast, nodes located towards the center of these clusters display 
increased n.s.i. cross-degree density, which is in general to be expected from the 
continuity of the underlying system. However, in that case we also note tendencies 
for decreased values of n.s.i. local cross-clustering coefficients. This observation is a 
result of the fact that especially those nodes with only one link to the opposite field 
show by definition the highest value of the n.s.i. local cross-clustering coefficient, 
C J V * = 1. With increasing n.s.i. cross-degree density this measure converges to a 
more reasonable estimate of a node’s tendency to cluster. 

One way to address this issue in the future would be to subtract the squared sum 
of weights w v of all neighbors of the considered node from the numerator in Eq. (10). 
Such procedure would, however, introduce a non-standard network measure whose 
properties should be assessed thoroughly in future research before applying it to 
climatic studies. To this end, we acknowledge that the concerned nodes do not 
play a crucial role for the propositions put forward in this section, since they (i) 
are ultimately dealt with by the assessment of n.s.i. cross-transitivity which weighs 
those corresponding nodes much lower than those with a high n.s.i. cross-degree 
density and (ii) only manifest in the very upper left corners of Fig. [6jA,B. In that 
case they do not contribute significantly to the observed relationship between n.s.i. 
cross-degree density and n.s.i. local cross-clustering coefficient and have no further 
impact on the qualitative statements put forward above. 

Following upon the quantitatively observed hierarchy, Fig. [7j allows for a visual 
inspection of some illustrative parts of the corresponding network structure. In 
particular, we display for two selected patches of nodes in the SST field that show 
high values of k 1 * with the HGT field (blue and orange shaded polygons in Fig. [?|4) 
their corresponding neighboring nodes in the HGT field as well as all links between 
those nodes (respective blue and orange scatter in Fig. [7)4). While ignoring very 
small clusters we find in total four (three) substantial mutually disconnected patches 
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of nodes in the HGT field that correlate with the respective ocean patches. Vice 
versa, by selecting the resulting two largest patches of nodes in the HGT field above 
both oceans (blue and orange shaded polygons in Fig. [?|3) we find that each of the 
patches only correlates with two disconnected patches of nodes in the SST field that 
are of a relevant spatial extent to have an effect on the estimation of C**. Thus, the 
resulting n.s.i. cross-clustering coefficient C 1 * for nodes in the ocean exceed C s * for 
nodes in the atmosphere as the ocean correlates with more mutually disconnected 
clusters of nodes than vice versa. 

Comparing the observed node patches in the HGT field with atmospheric pat¬ 
terns of large-scale variability patterns (Handorf and Dethloff 2012), we relate the 
two atmospheric clusters in the HGT field that are located above the Atlantic (blue 
scatter in Fig. [7]A) with the NAO. Correspondingly, the three patches located above 
the Pacific (orange scatter in Fig. [f}4) coincide well with the spatial signature of 
the PNA pattern. Taking into account past studies that applied lagged correlation 
analysis we note that on the time scales considered in this study the atmosphere 
serves as a driving force of the ocean along the spatial domain that is of interest 
here (e.g. Czaja and Frankignoul 1999; Frankignoul and Sennchael 2007; Gastineau 
and Frankignoul 2015). Thus, the hierarchical network structure might on the one 
hand be a result of the aforementioned atmospheric forcing to the ocean. On the 
other hand, with reference to Fig. [7| the framework of coupled climate networks 
and the methodology put forward in this work serve to resolve the corresponding 
induced correlation structure between the two climatic subsystems in a spatially 
explicit way, such that it enables to specifically detect forcing and forced areas in 
atmosphere and ocean, respectively. 

Choosing different HGT layers up to 200 mbar yield similar results (not shown). 
This aligns well with previous results by Czaja and Frankignoul (1999) and Frankig¬ 
noul and Sennchael (2007) who observed comparable spatial statistical patterns at 
each tropospheric level. Thus, the observed hierarchy seems to be a generic property 
of the troposphere. Following the above lines of thought, future work should investi¬ 
gate coupled climate networks constructed from lagged cross-correlations to investi¬ 
gate whether the observed structures are indeed a result of short-term atmosphere- 
to-ocean forcing. Such procedure would, however, require the derivation of novel 
directed interacting network measures, which in turn would provide a valuable ex¬ 
tension to the framework of climate network analysis. 


3.3 Global measures 

So far we have focused our study on three atmospheric layers, namely the 50 mbar, 
100 mbar and 500 mbar HGT field. Specifically for the latter case, we have carried 
out a further in-depth analysis of the observed hierarchical structures by means of 
assessing the power-law relationship between k 1 * and C 1 * as well as investigating the 
distinct spatial distribution of nodes and links in ocean and atmosphere that obey 
the observed hierarchical organization (Fig. [Tj) . To show that these structures are 
(i) not only present for the 500 mbar field and (ii) their observations are robust with 
respect to the choice of fink densities we investigate global network characteristics 
that provide a macroscopic description of the observed network structures. 
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Specifically, we study the n.s.i. cross-transitivity T* computed over nodes in 
the SST field and T*. computed over nodes in each of the HGT fields according 
to Eq. (16) together with the n.s.i. global cross-clustering coefficients C* si and C* s , 
respectively (Eq. ®). Note again that the latter are defined as the weighted means 
of their local counterparts that are presented in Fig. [4j where nodes with no links 
to the opposite field are weighted in the same fashion as those with adjacent cross¬ 
links. In contrast, the n.s.i. cross-transitivity assigns nodes a weight corresponding 
to their n.s.i. cross-degree (which for the case of no adjacent cross-links takes a value 
of zero) and, thus, excludes them from the averaging. 

The corresponding results are summarized in Fig. [8} We find that both 1~* and 
C* si show their maximum values at around 10 km altitude (250 mbar) (Figs. [8}4 
and [8^3). For the same quantities, distinct minima at 850 mbar (1.4 km) coincide 
with the transition from the atmospheric boundary layer to the lower troposphere 
as also found in Donges et al. (2011). For all layers above 100 mbar, T* remains 
almost constant at low values. Hence, T* and C* sl seem to naturally discriminate 
between three different atmospheric layers: below 850 mbar (atmospheric boundary 
layer), between 850 mbar and 100 mbar (free troposphere) and above 100 mbar 
(lower stratosphere). 

For the global measures computed over all nodes in the HGT field, we find that 
the n.s.i. cross-transitivity T* s shows almost constant values for all layers below 200 
mbar and, hence, again separates well the dynamics within the troposphere from 
that inside the stratosphere (Fig. |8p). For all layers above 200 mbar 7~* becomes 
almost independent of the cross-link density p s i that is fixed when constructing the 
network. The same property also holds for the n.s.i. global cross-clustering coefficient 
C* s computed over all nodes in the different HGT fields (Fig. |8p). 

In agreement with the local measures discussed in Sec. |3.2| we find that n.s.i. cross¬ 
transitivity and n.s.i. global cross-clustering coefficients are in most cases larger for 
nodes in the HGT fields than for nodes in the SST field (compare Fig. [8]4 with 
Fig. |8p and Fig. [8j3 with Fig. [8jD) . As for the n.s.i. local cross-clustering coefficients 
this indicates again the hierarchical network structure, i.e., a tendency for nodes 
in the HGT field to form triangular structures with nodes in the SST field, that 
is present across all atmospheric layers ranging from the troposphere to the lower 
stratosphere. The detailed structure of this hierarchy, however, seems to vary with 
the different atmospheric layers under study. 

This observation further holds not only for the case of p si = 0.005 that was used 
in the previous sections but also for larger values that are chosen from a reasonable 
range (Fig. |8j). Thus, we consider our results to be sufficiently robust with respect 
to the choice of the networks’ link densities. 

In general, we observe that the quantitative and qualitative properties of the 
n.s.i. cross-transitivity and n.s.i. global cross-clustering coefficients vary with the 
different atmospheric layers. Hence, these global characteristics may serve to inter¬ 
compare and distinguish between different correlation structures in a coupled climate 
network. An in-depth analysis of the mechanisms that cause the occurrence of this 
behavior in our specific application remains as a subject of future research. 
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4 Conclusions & Outlook 


We have carried out a detailed analysis of the correlation structure between at¬ 
mospheric and ocean dynamics in the Northern Hemisphere extratropics from the 
viewpoint of coupled climate networks. Comparison between the n.s.i. cross-degree 
density (measuring the weighted share of significant correlations between grid points 
in different layers) and the leading mode of the maximum covariance analysis (MCA) 
reveals an expected high congruence between both methods for the considered data 
sets. However, coupled network analysis, and particularly the investigation of higher- 
order network parameters, allows us to further disentangle the underlying correlation 
structure. The (average) n.s.i. cross-degree density in combination with the (aver¬ 
age) n.s.i. local cross-clustering coefficient provides additional insights on areas in 
the ocean and the atmosphere that show strong mutual correlations as well as local¬ 
ized versus delocalized correlation structures with the respective opposite field. In 
the SST field nodes tend to correlate with multiple mutually unconnected groups of 
similar nodes within the respective HGT fields. From investigating the interdepen¬ 
dency between n.s.i. cross-degree density and n.s.i. local cross-clustering coefficient, 
we have found that the correlations between the ocean and the atmosphere exhibit a 
hierarchical structure in the sense of a power-law relationship between both proper¬ 
ties. A visual inspection of the coupled climate network for the case of the 500 mbar 
HGT field reveals that the observed structure could be related with a forcing of the 
ocean by the dominant atmospheric patterns above the Atlantic and the Pacific. 
Ultimately, global network characteristics further support the results obtained from 
their local correspondents by showing that the observed structure is valid for large 
parts of the atmosphere ranging from the troposphere to the lower stratosphere. 

In order to discriminate between the internal variability of the fields under study 
and possible influences of an external forcing, future work should analyze ensem¬ 
ble simulations of general circulation models to rule out common driver effects or 
assess the likelihood of their influence on the observed structures. In order to in¬ 
vestigate the influence of spatio-temporal auto-correlation on the outcome of the 
present analysis the network could be alternatively constructed by estimating pair¬ 
wise thresholds from surrogate data as proposed by Palus et al. 2011. This approach 
would, however, break the comparability of the network approach with that of max¬ 
imum covariance analysis, such that a different way of validating and comparing 
the results must be found. Comparability could be achieved by assessing synthetic 
model data, e.g., created from an auto-regressive process based on principal compo¬ 
nents of the data sets under study, and the application of both, MCA and network 
analysis. In addition to probable influences of auto-correlation such a process would 
allow to assess the influence of different time scales in ocean and atmosphere on the 
involved network characteristics. 

Besides all possible future lines of work with respect to the climatic side of this 
work, from a network-theoretic point of view it is worthwhile to construct climate 
networks using more advanced causal estimators (e.g. Runge et al. 2012; Runge et 
al. 2014; Runge et al. 2015) to disentangle direct from indirect or externally and 
internally forced correlations. 

To this end, our analysis has only been performed for the pairwise correlation 
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between one atmospheric layer and the ocean. Future studies should further ex¬ 
plore the possibility to refine the proposed methods to also quantify interactions in 
a climate network existing of more than two subnetworks. Specifically, when study¬ 
ing coupled climate networks in the Northern Hemisphere, one should also consider 
Arctic sea ice as an additional observable in the network construction. Its dynamics 
has already been discovered as an influencing factor on atmospheric teleconnections 
and the dynamics of land snow cover in the Northern Hemisphere (Handorf et al. 


2015). The study of coupled climate networks can help here to further disentangle 


and quantify possible changes in correlations between ocean and atmosphere over 
the course of the past decades that may have been induced by processes related to 
the Arctic amplification (Serreze and Francis 2006). Moreover, it is of great interest 
to apply our methods not only to coupled networks composed of different climatic 
fields (as presented in this work), but also to networks constructed from just one 
single climatic field that divides into dynamically distinct areas (Hlinka et al. 2014) 
or communities (Steinhaeuser et al. 2011; Tsonis et al. 2010). The framework pre¬ 
sented in this work could then be utilized to study and quantify correlations between 
these detected or defined regions on or parallel to the Earth’s surface. This would 
allow for a detailed investigation of correlation structures between different climatic 
subsystems such as, for example, the Indian Summer Monsoon and the El Nino 
Southern Oscillation. 

Finally, it remains to remark that the weighted network measures presented in 
this work provide a general framework which can be applied to quantify interdepen¬ 
dencies in complex networks representing subjects of study taken from many other 
fields beyond climatology. 
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Figure 1: Cross-threshold T si between the subnetwork constructed from the SST field 
and all 18 isobaric surfaces of HGT in winter for different standard (unweighted) 
cross-link densities. 
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Figure 2: Leading coupled patterns obtained from MCA between the SST field and 
three layers of geopotential height at 50 mbar (A and B), 100 mbar (C and D) and 
500 mbar (E and F) in winter (DJF). The left column (A, C and E) displays the 
component in the SST and the right column (B, D and F) that in the respective 
HGT field. All spatial patterns are shown as regression maps obtained by regressing 
SST anomalies and geopotential height anomalies onto associated time series for 
the geopotential height field for the respective MCA mode. Statistically significant 
areas at the 95% confidence level based on a two-tailed Student’s t-test are shown 
as black contours (for the SST maps) and grey contours (for the geopotential height 
maps). Dashed lines indicate the climatic mean geopotential height fields. 
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Figure 3: N.s.i. cross-degree density for coupled climate networks constructed from 
the SST field and three layers of geopotential height at 50 mbar (A and B), 100 mbar 
(C and D) and 500 mbar (E and F) for winter months (DJF). The left column (A, 
C and E) displays the n.s.i. cross degree density n l * for links pointing from the SST 
into the HGT subnetwork while the right column (B, D and F) displays the n.s.i. 
cross-degree density k s * for links pointing from the HGT into the SST subnetwork. 
Only nodes with n l * > 0 and k s * > 0 are shown. 


27 























































150°W 180° 150°E 



150°W 180° 150°E 

120 ° 

90°W 

60° 

30°W 0° 30°E 

150°W 180° 150°E 

120 ° 

90°W 

60° 

30°W 0° 30°E 






150°W 180° 150°E 


30°W 0° 30°E 

150°W 180° 150°E 


30°W 0° 30°E 



120 


90 


60 


150°W 180° 150°E 


30°W 0° 30°E 


1.0 

0.8 


0.6 


0.4 




0.2 


0.0 


1.0 

0.8 


0.6 


0.4 


* 

CO 




0.2 


0.0 


1.0 

0.8 


0.6 


0.4 


* 

CO 




0.2 


0.0 


Figure 4: As in Fig. [ 3 ] for the n.s.i. local cross-clustering coefficients Q* and C**. 
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Figure 5: Schematic explanation of the observed quantitative differences in the n.s.i. 
local cross-clustering coefficients for nodes in the SST and HGT fields. Nodes in 
the ocean (box 3) tend to connect with statistically dissimilar and thus unconnected 
clusters of nodes in the atmosphere (such as nodes in box 1 and 2). Hence, the 
n.s.i. local cross-clustering coefficient Q* only takes low values. In contrast, nodes 
in the atmosphere, e.g. from box 1, likely connect with clusters in the SST field, 
such as nodes exclusively in box 3. This results in a high n.s.i. cross local-clustering 
coefficient C** for nodes in the atmosphere. 
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Figure 6: N.s.i. local cross-clustering coefficients C**(/%**) for nodes in the SST field 
(A) and CI*(k, s *) for nodes in the 500 mbar HGT field (B) as functions of the re¬ 
spective n.s.i. cross-degree densities. The dashed fine in (A) indicates the relation¬ 
ship C l * ~ (k**) “ (here with a = 0.94) expected for traditional network measures 
C v (k v ) in the case of hierarchical network structures (Ravasz et al. 2002; Ravasz and 
Barabasi 2003). 
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Figure 7: Visualization of a selection of nodes that are relevant for the observed 
hierarchical network structure. (A) Two clusters of nodes in the SST field (blue and 
orange shaded polygons) that show positive values of k 1 * with the 500 mbar HGT 
field (compare Fig. [3^). Correspondingly, coloured scatter points denote nodes in 
the HGT field that the considered SST nodes are connected with. All links that 
mutually connect the resulting HGT nodes are displayed as well. (B) The same for 
the two largest patches of nodes in the 500 mbar HGT field that where detected 
in (A). Coloured scatter points now indicate all nodes in the SST field that are 
connected with these patches. 
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Figure 8: Global coupled network measures computed for all 18 coupled climate 
networks: (A) N.s.i. cross-transitivity and (B) n.s.i. global cross-clustering coeffi¬ 
cient (B) taken over all nodes in the SST field. (C) and (D) display the respective 
measures computed over all nodes in the HGT field. To demonstrate the robustness 
and consistency of the results, we construct the networks for different choices of 
(unweighted) cross-link density p si and internal link density Pi = p s = ‘^Psi- 
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Table 1: Air pressure p, and associated mean geopotential height Z, as well as 
the internal threshold T, ( p vl = 0.01) corresponding to an internal link density of 
pa = 0.01 for each isobaric surface i. 

Layer Air pressure Geopotential Threshold 
i Pi [mbar] height Zi [km] Ti{pa = 0.01) 


0 

10 

30.9 

0.9919 

1 

20 

26.3 

0.9936 

2 

30 

23.7 

0.9932 

3 

50 

20.5 

0.9876 

4 

70 

18.4 

0.9781 

5 

100 

16.2 

0.9621 

6 

150 

13.7 

0.9263 

7 

200 

11.8 

0.9166 

8 

250 

10.4 

0.8982 

9 

300 

9.2 

0.8894 

10 

400 

7.2 

0.8895 

11 

500 

5.6 

0.8958 

12 

600 

4.2 

0.9036 

13 

700 

3.0 

0.9119 

14 

775 

2.2 

0.9171 

15 

850 

1.4 

0.9205 

16 

925 

0.8 

0.9215 

17 

1000 

0.1 

0.9197 
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